
* Full employment labor market CGE model

* Produces results for central calibration carbon tax with lump-sum rebates

* Define sets 

sets
        t "time periods" /t1*t51/
        tp(t) /t1*t51/
        tpi(t) /t1/
        tfirst(t) "First time period"
        tlast(t) "Final time period"
        ind /i1*i22/
        ie(ind) "Energy Industries" /i1*i2,i5*i7/
        im(ind) "Material Industries" /i3,i4,i8*i22/
        oilgas(ind) "Oil&Gas Extraction" /i1/
        coal(ind) "Coal Mining" /i2/
        industrial(ind) "Industrial" /i1*i4,i8*i12/
        commercial(ind) "Commercial" /i13*i14,i16*i21/
        gov(ind) "Government" /i22/
        transportation(ind) "Transportation" /i15/
        elec(ind) "Elecric Power" /i5/ 
        ngd(ind) "Nat Gas Dist" /i6/
        refine(ind) "Refine" /i7/
        refine_n(ind) /i1*i6,i8*i22/
        primary(ind) "Primary Energy" /i1*i2/ 
        secondary(ind) "Secondary Energy" /i6*i7/ 
        durable(ind) "Durable Manuf" /i11/ 
        nondurable(ind)  "Non-durable Manuf" /i12/
        other(ind) /i3,i4,i8,i9,i10,i16/
        trade(ind) "Trade"  /i13,i14/
        it "iterations" /it1*it20/
        out /out1*out12/
        services(ind) /i17*i21/
        i /i1*i22,Gov,All/    
;

tfirst(t) = yes$(ord(t) eq 1);
tlast(t) = yes$(ord(t) eq card(t));
alias (ind,ind2)
alias (ind,ind3)
alias(ie,ie2);
alias(im,im2);
alias(industrial,industrial2)
alias(commercial,commercial2)
alias(refine_n,refine_n2)



* DEFINE DEEP PARAMETERS
parameters
         beta Time preference" /.996/
         chi "Frisch elasticity of labor supply" /1/
         psi "Disutility of work"
         psi_f "Disutility of work - foreign"
         mu_e(ind,ind2) "Emissions factor"
         tau_l0 "Tax on labor income" /0.31/
         tau_p0 "Payroll tax" /0.06/
         sigma /1/
         sigmac "Elast. of subs. between cons. goods" /0.75/
         sig_e(ind) "Elast. of subs. between energy goods"
         sig_eg "Elast. of subs. between energy goods - gov"
         sig_m(ind) "Elast. of subs. between non-energy goods"
         sig_mg "Elast. of subs. between non-energy goods - gov"
         sig_v(ind) "Elast. of subs. between energy and non-energy goods"
         sig_vg "Elast. of subs. between energy and non-energy goods - gov"
         sig_y(ind) "Elast. of subs. between labor and inter. inputs"
         sig_yg "Elast. of subs. between labor and inter. inputs - gov"
         rho(ind)
         rhog

         sig_df(ind)
         sig_c_f(ind)
         sig_c(ind)


;

sig_y(ind) = 0.5;
sig_e(ind) = 0.4;
sig_m(ind) = 0.6;
sig_v(ind) = 0.7;

sig_v('i5')  = 0.4;
sig_v('i6')  = 0.4;
sig_v('i7')  = 0.4;
sig_v('i8')  = 0.4;

sig_eg = 0.2;
sig_mg = 0.2;
sig_vg = 0.2;
sig_yg = 0.2;

sig_e(elec) = 1.5;

rho(ind) = (sig_y(ind)-1)/sig_y(ind);
rhog = (sig_yg-1)/sig_yg;

sig_df(ind) = 1.5;
sig_c(ind) = 1.5;
sig_c_f(ind) = 1.5;




$INCLUDE FECGE_JAERE_CAL.gms


* DEFINE STEADY STATE CARBON PRICE
parameter
p_e0 /0/
tar_i2(ind)
;

tar_i2(refine) = f_y0(refine)*sum(oilgas,mu_e(oilgas,refine));

$INCLUDE FECGE_JAERE_EQN.gms



* REPLICATE BASE CASE WITH P_E = 0;
options iterlim = 0;
delta_tau.fx = 1;

solve fecge_ss_ct using mcp ;


* SOLVE STEADY STATE WITH P_E>0
options iterlim = 10000000;

* Define steady state output parameters

parameter
l_d
l_d_perc
tlump_d
net_gross
crev_ss_ct
perc_e
perc_l
utot_ss_ct
tbal_db
mc_db
bc_db
lmc_db
mc_f_db
bc_f_db
lmc_f_db
yg_db
yg_f_db
perc_imports
perc_exports
lg_d
lg_d_perc
perc_l_new
p_e
perc_w
perc_income
p_f_ss_perc
;



loop(it,p_e0 = p_e0+5;



solve fecge_ss_ct using mcp  ;


* SAVE STEADY STATE OUTPUT


p_e(it) = p_e0;
l_d(ind)= l_ss_ct.l(ind)-l_ss(ind);
l_d_perc(ind,it)= (l_ss_ct.l(ind)-l_ss(ind))/l_ss(ind);

lg_d = lg_ss_ct.l - lg_ss;
lg_d_perc(it) = lg_d/lg_ss;

tlump_d(it) = tlump_ss_ct.l - tlump_ss;
crev_ss_ct(it) = ( sum(ie,p_e0*sum(ind,mu_e(ie,ind)*(iod_ss_ct.l(ie,ind)+iof_ss_ct.l(ie,ind)))) + sum(ind,tar_i2(ind)*p_e0*(sum(ind2,iof_ss_ct.l(ind,ind2))+cf_ss_ct.l(ind))) )/p_cbar_ss_ct.l;
net_gross(it) = tlump_d(it)/crev_ss_ct(it);

perc_e(it) = (etot_ss_ct.l - etot0)/etot0;
perc_l(it) = (sum(ind,l_ss_ct.l(ind))+lg_ss_ct.l-(sum(ind,l_ss(ind))+lg_ss))/(sum(ind,l_ss(ind))+lg_ss);
utot_ss_ct(it) = 1 - 0.95*(1+perc_l(it));
perc_w(it) = (w_ss_ct.l/p_cbar_ss_ct.l)/(w_ss/p_cbar_ss) - 1;
perc_income(it) = (w_ss_ct.l*ls_ss_ct.l/p_cbar_ss_ct.l)/(ltot_ss*w_ss/p_cbar_ss)-1;

tbal_db(it) = sum(ind,p_ss_ct.l(ind)*exports_ss_ct.l(ind) - (p_f_ss_ct.l(ind)/exch_ss_ct.l)*imports_ss_ct.l(ind)) ;
mc_db(ind,it) = y_ss_ct.l(ind) - cd_ss_ct.l(ind) - gd_ss_ct.l(ind) -  exports_ss_ct.l(ind) - sum(ind2,iod_ss_ct.l(ind,ind2)) ;
bc_db(it) = (1-tau_l0)*ls_ss_ct.l*w_ss_ct.l + p_cbar_ss_ct.l*tlump_ss_ct.l  + sum(ind,prof_ss_ct.l(ind))  - p_cbar_ss_ct.l*cbar_ss_ct.l;
lmc_db(it) = ls_ss_ct.l - sum(ind,l_ss_ct.l(ind))-lg_ss_ct.l;
mc_f_db(ind,it) = y_f_ss_ct.l(ind) - cd_f_ss_ct.l(ind) - gd_f_ss_ct.l(ind) -  imports_ss_ct.l(ind) - sum(ind2,iod_f_ss_ct.l(ind,ind2)) ;
bc_f_db(it) = (1-tau_l0)*ls_f_ss_ct.l*w_f_ss_ct.l + p_cbar_f_ss_ct.l*tlump_f_ss_ct.l  + sum(ind,prof_f_ss_ct.l(ind))  - p_cbar_f_ss_ct.l*cbar_f_ss_ct.l;
lmc_f_db(it) = ls_f_ss_ct.l - sum(ind,l_f_ss_ct.l(ind))-lg_f_ss_ct.l;
yg_db(it) =  gammayg0*( alphapg*lg_ss_ct.l**rhog + (1-alphapg)*vvg_ss_ct.l**rhog)**(1/rhog) - yg0;
yg_f_db(it) =  gammayg0_f*( alphapg*lg_f_ss_ct.l**rhog + (1-alphapg)*vvg_f_ss_ct.l**rhog)**(1/rhog) - yg0_f;

loop(ind,
if(exports0(ind)>0,
perc_exports(ind,it) = (exports_ss_ct.l(ind)-exports0(ind))/exports0(ind);
);
if(imports0(ind)>0,
perc_imports(ind,it) = (imports_ss_ct.l(ind)-imports0(ind))/imports0(ind);
);
);

p_f_ss_perc(ind,it) = p_f_ss_ct.l(ind)/exch_ss_ct.l;
);


* DISPLAY KEY RESULTS
options decimals = 5;


display
l_d_perc
lg_d_perc
net_gross
perc_e
perc_l
utot_ss_ct
tbal_db
mc_db
bc_db
lmc_db
mc_f_db
bc_f_db
lmc_f_db
yg_db
yg_f_db
perc_exports
perc_imports
p_f_ss_perc
exch_ss_ct.l
p_ss_ct.l
w_ss_ct.l
p_f_ss_ct.l
w_f_ss_ct.l

;


parameter
figure_FECGE_LUMP(*,it)
table_FECGE_LUMP_p_e40(*,*)
;

figure_FECGE_LUMP('out1',it) = -1*perc_e(it);
figure_FECGE_LUMP('out2',it) = p_e(it);
figure_FECGE_LUMP('out3',it) = utot_ss_ct(it)-.05;
figure_FECGE_LUMP('out4',it) = perc_l(it);
figure_FECGE_LUMP('out5',it) = perc_w(it);
figure_FECGE_LUMP('out6',it) = perc_income(it);


table_FECGE_LUMP_p_e40(ind,'out1') = l_d_perc(ind,"it8");
table_FECGE_LUMP_p_e40('Gov','out1') = lg_d_perc("it8");
table_FECGE_LUMP_p_e40('All','out1') = perc_l("it8");




execute_unload "table_JAERE_FECGE_LUMP_p_e40.gdx"     table_FECGE_LUMP_p_e40
execute 'gdxxrw.exe table_JAERE_FECGE_LUMP_p_e40.gdx  o=table_JAERE_FECGE_LUMP_p_e40.xls SQ=N par=table_FECGE_LUMP_p_e40 rng=table_FECGE_LUMP_p_e40!'

execute_unload "figure_JAERE_FECGE_LUMP.gdx"     figure_FECGE_LUMP
execute 'gdxxrw.exe figure_JAERE_FECGE_LUMP.gdx  o=figure_JAERE_FECGE_LUMP.xls SQ=N par=figure_FECGE_LUMP rng=figure_FECGE_LUMP!'

display
figure_FECGE_LUMP
table_FECGE_LUMP_p_e40
;

